Design of viscoelastic coatings to reduce turbulent friction drag

ABSTRACT

A method is provided to select appropriate material properties for turbulent friction drag reduction, given a specific body configuration and freestream velocity. The method is based on a mathematical description of the balance of energy at the interface between the viscoelastic surface and the moving fluid, and permits determination of the interaction of turbulent boundary layer fluctuations with a viscoelastic layer by solving two subtasks—i.e., a hydrodynamic problem and an elasticity problem, which are coupled by absorption and compliancy coefficients. Displacement, velocity, and energy transfer boundary conditions on a viscoelastic surface are determined, and a Reynolds stress type turbulence model is modified to account for redistribution of turbulent energy in the near-wall region of the boundary layer. The invention permits drag reduction by a coating with specified density, thickness, and complex shear modulus to be predicted for a given body geometry and freestream velocity. For practical applications, viscoelastic coatings may be combined with additional structure, including underlying wedges to minimize edge effects for coatings of finite length, and surface riblets, for stabilization of longitudinal vortices.

This application claims the benefit of No. 60/138,023, filed Jun. 8,1999.

BACKGROUND OF INVENTION

Since M. O. Kramer reported successful experimental results in 1957,there have been repeated attempts to reduce frictional drag in turbulentfluid flow over a surface by applying a passive compliant coating.Experimental results in this area have been mixed. Most investigatorshave reported a drag increase, while only a few have claimed dragreduction for turbulent flow. A number of theoretical studies havecharacterized the stability of the laminar boundary layer over adeforming surface and other studies have characterized the reaction of acoating to a fluctuating load. However, no rigorous analytical techniquehas been previously reported that has been used to successfully design adrag-reducing coating for turbulent flow.

In the past, passive coatings were tested without specification and fullcharacterization of critical physical parameters, such as the frequencydependent complex shear modulus, density, and thickness. In order toachieve and ensure drag reduction with a viscoelastic coating, amethodology is required for selecting appropriate material propertiesand for estimating anticipated drag reduction as a function ofconfiguration and velocity.

Relevant background information for associated technical topics isavailable in the literature, and may be useful due to the technicalcomplexity of this invention. A classical discussion of boundary layertheory, including formulation of Navier-Stokes and turbulent boundarylayer equations, is provided in Boundary-Layer Theory, by Dr. HermannSchlichting, published by McGraw Hill, New York, seventh edition, 1979.A discussion of structures and scales in turbulent flows can be found inTurbulence, 1975, McGraw Hill, written by J. O. Hinze, and in “CoherentMotions in the Turbulent Boundary Layer,” in Annual Review of FluidMechanics, 1991, volume 23, pp. 601-39, written by Steven K. Robinson.Background on Reynolds stress types of turbulence models is found in thechapter, “Turbulent Flows: Model Equations and Solution Methodology,”written by Tom Gatski, and included in the Handbook of ComputationalFluid Mechanics, published by Academic Press in 1996. Equations in fluidand solid mechanics are often expressed in indicial, or tensor,notation, for compactness. Chapter 2 in the text A First Course inContinuum Mechanics, by Y. C. Fung, Prentice-Hall, Inc., EnglewoodCliffs, N.J., 1977, provides a brief introduction into tensor notationfor mechanics equations. An introduction to finite difference methods,which are used to solve the system of momentum and continuity equationsfor a turbulent fluid, is provided in the text, Computational FluidDynamics for Engineers, written by Klaus Hoffman, and published in 1989by the Engineering Education System i Austin in Austin, Texas.Descriptions of measured and mathematically modeled physical propertiesof polymers are found in the text, Viscoelastic Properties of Polymersby J. D. Ferry, Wiley, New York, 1980, 3^(rd) edition. The article,“Loss Factor Height and Width Limits for Polymer Relaxation,” by BruceHartmann, Gilbert Lee, and John Lee, in the Journal of the AcousticalSociety of America Vol. 95, No. 1, January 1994, discusses mathematicalcharacterization of shear moduli for real viscoelastic, polymericmaterials, including those approximated by the Havriliak-Negamiapproach.

Recently in the international literature (K. S. Choi, X. Yang, B. R.Clayton, E. J. Glover, M. Atlar, B. N. Semonev, and V. M. Kulik,“Turbulent Drag Reduction Using Compliant Surfaces,” Proceedings of theRoyal Society of London, A (1997) 453, pp. 2229-2240). Choi et al.reported experimental measurements of up to 7% turbulent friction dragreduction for an axisymmetric body coated with a viscoelastic material.These experiments were performed in the United Kingdom, using coatingsdesigned and fabricated in Russia at the Institute of Thermophysics,Russian Academy of Sciences, Novosibirsk, by a team headed by B .N.Semenov. The basic design approach was outlined in “On Conditions ofModelling and Choice of Viscoelastic Coatings for Drag Reduction,” inRecent Developments in Turbulence Management, K. S. Choi, ed., 1991, pp.241-262, Dordrecht, Kluwer Publishers. The Novosibirsk design approachis semi-empirical in nature, and does not take into account the fullcharacterization of the complex shear modulus of the viscoelasticmaterial, namely, the relaxation time of the material. The Novosibirskdesign approach does take into account frequency-dependent materialproperties. Furthermore, the Novosibirsk concept is valid only for amembrane-type coating, such as a film which coats a foam-rubbersaturated with water or glycerine, and where only normal fluctuations ofthe surface are considered.

The structure of coatings intended for drag reduction has been addressedin the international literature, starting with the 1938 patent No.669-897, “An Apparatus for the Reduction of Friction Drag,” issued inGermany to Max O. Kramer. Kramer later received a patent in 1964, U.S.Pat. No. 3,161,385, and in 1971, U.S. Pat. No. 3,585,953 for coatings toextend laminar flow in a boundary layer. Soviet inventor's certificates,such as “A Damping Covering,” USSR patent 1413286, Publication20.01.1974, Bulletin of the Inventions 14, by V. V. Babenko, L. F.Kozlov, and S. V. Pershin, “An Adjustable Damping Covering,” USSR patent1597866, Publication 15.03.1978, Bulletin of the Inventions 110, by V.V. Babenko, L. F. Kozlov, and V. I. Korobov, and “A Damping Covering forSolid Bodies,” USSR patent 1802672, Publication 07.02.1981, Bulletin ofthe Inventions 15, by V. V. Babenko and N. F. Yurchenko, have alsodescribed the structure of drag-reducing coatings comprised ofviscoelastic materials. These inventor's certificates identified thethree-dimensional structure within a drag-reducing coating, but do notaddress the methodology for choosing appropriate parameters of theviscoelastic materials to be used in the manufacture of such coatings.Structural features include multiple layers of materials, longitudinal,rib-like inclusions of elastic, viscoelastic, or fluid materials, andheated elements. Viscoelastic coatings may be combined with other formsof structure, such as longitudinal riblets molded on or within thesurface of the coating. As described in the international literature inpublications such as “Secondary Flow Induced by Riblets,” written by D.B. Goldstein and T. C. Tuan, and published in the Journal of FluidMechanics, volume 363, May 25, 1998, pp. 115-152, two-dimensional, rigidriblets alone have been shown experimentally to reduce surface frictiondrag up to about 10%.

BRIEF SUMMARY OF THE INVENTION

The present invention enables the design of a passive viscoelasticcoating for the reduction of turbulent friction drag. Coatings withmaterial properties designed using the methodology described in thisinvention have reduced friction drag by greater than 10%. Themethodology of the present invention permits, as a first object of theinvention, the specification of the frequency dependent complex shearmodulus, the density, and the thickness of an isotropic viscoelasticmaterial which will reduce turbulent friction drag relative to specificflow conditions over a rigid surface. Quantitative levels of dragreduction can be estimated. Mathematical detail is provided for thecases of turbulent flow over a rigid flat plate as well as aviscoelastic flat plate, where the invention accounts for both normaland longitudinal oscillations of the surface. A second object of theinvention is the specification of material properties for a coatingcomposed of multiple layers of isotropic viscoelastic materials. A thirdobject of the invention is the specification of material properties fora coating composed of an anisotropic material. A fourth object of theinvention is the minimization of edge effects for coatings of finitelength. A fifth object of the invention is the stabilization oflongitudinal vortices through combination of viscoelastic coating designwith additional structure, such as riblets.

The methodology used herein to describe the interaction of a turbulentboundary layer (TBL) with a viscoelastic (VE) layer involves twotasks, 1) a fluids task, involving the calculation of turbulent boundarylayer parameters, given boundary conditions for a rigid, elastic, orviscoelastic surface (herein referred to as the TBL problem), and 2) amaterials task, involving the calculation of the response of aviscoelastic or elastic surface to a periodic forcing function whichapproximates the loading of the turbulent boundary layer. The inventionfocuses upon cation of amplitudes of layer. The invention focuses uponcalculation of amplitudes of surface oscillations and velocities, and ofthe energy flux for a viscoelastic coating (hereinafter referred to asthe VE problem). These two tasks are coupled by coefficients related tosurface boundary conditions of energy absorption and surface oscillationamplitudes (hereinafter referred to as dynamic and kinematic boundaryconditions, respectively). The TBL problem is first solved for a rigidsurface, thus providing necessary input to describe the forcing functionon the surface, and also providing baseline calculations of frictiondrag, for comparison. The VE problem is solved next, given a periodicforcing function that approximates the shear and pressure pulsations ofa given boundary layer. Initial choices for material parameters arebased on theoretical and empirical guidelines. Optimal materialparameters are chosen, following a series of iterations, such that thefollowing two criteria are met:

-   -   1) the energy flux into the viscoelastic coating is maximum, and    -   2) the amplitudes of surface oscillation are less than the        thickness of the viscous sublayer of the turbulent boundary        layer over the coating.        If the amplitude of oscillations exceeds the thickness of the        viscous sublayer of the turbulent flow, then the oscillations        effectively increase the roughness of the surface, thus leading        to an increase in friction drag. Furthermore, as the phase speed        of disturbances in the boundary layer exceeds the shear wave        speed in the material, a resonant interaction with large        amplitude waves occurs. These conditions are to be avoided.        Moderate energy flux into the material, however, where energy is        transformed into internal shear waves and eventually dissipated        as heat, leads to a qualitative and quantitative change in the        turbulent energy balance, with a consequent reduction in        friction drag. Optimal physical properties of a viscoelastic        material will vary with freestream velocity, position along a        body, pressure gradient, and any other factors which influence        the development of the boundary layer and the characteristics of        local turbulent fluctuations.

By solving the TBL equations, the turbulent friction drag over aviscoelastic, elastic, or rigid surface can be quantitatively evaluated.In the case of a viscoelastic surface, where energy is absorbed andsurface oscillations are nonzero, both dynamic and kinematic boundaryconditions are specified. These boundary conditions are derived directlyfrom the solution of the VE equations for energy flux and surfaceoscillation amplitudes, and then transferred into a dissipation boundarycondition and Reynolds stress boundary conditions for solution of theTBL equations. Vertical oscillations influence the effective roughnessof the surface, and the root-mean square (rms) value of the verticaloscillation amplitude is classified as the dynamic roughness. If theoscillation amplitudes are lower than the viscous sublayer thickness, itis appropriate to estimate Reynolds stresses as zero. The equations fora turbulent boundary layer describe turbulent diffusion as a gradientapproximation, which accommodates the dynamic boundary condition, andnear-wall functions are introduced to describe for different surfacesthe redistribution of turbulent energy in the near-wall region.

BRIEF DESCRIPTION OF THE DRAWINGS

The present invention will become more fully understood from thedetailed description given below and the accompanying drawings, whichare given by way of illustration only and thus are not limitative of thepresent invention, wherein:

FIG. 1a is a schematic of the interaction of a passive viscoelasticsurface with a turbulent boundary layer, where a low-amplitude travelingwave develops on the surface, FIG. 1b is an enlargement of the regionwithin the circle shown in FIG. 1a, and indicates displacementamplitudes, and

FIG. 2 is a flowchart of the basic methodology required to a) choose aviscoelastic material that reduces turbulent friction drag under givenflow conditions, and b) to quantify the value of drag reduction relativeto a rigid surface, and

FIG. 3 is a detail of a structure introduced at the intersection betweena viscoelastic coating and a rigid surface in order to minimize localoscillation amplitude and edge effects.

DETAILED DESCRIPTION

The present invention identifies physical and geometric parameters of aviscoelastic coating that reduces turbulent friction drag under givenflow conditions. Furthermore, the invention permits evaluation of theanticipated drag reduction effectiveness of a given material with knownphysical properties for a given body configuration and set of flowconditions. The methodology has been applied principally to thecharacterization of coatings for turbulent flow over flat plates andbodies of revolution, and can also be applied to more complex geometrieshaving curvature and nonzero pressure gradients.

A fluid boundary layer is the very thin layer of fluid adjacent to asurface over which fluid is flowing. It is the region where frictionalforces play a major role, and is where the flow adjusts from conditionsat the surface to conditions in the freestream of the flow. The outeredge of the boundary layer is traditionally defined as that locationwhere the ratio, β, of the mean velocity, U, to the freestream velocity,U₂₈ , is a constant which is approximately equal to 1: $\begin{matrix}{{\frac{U}{U_{\infty}}}_{\frac{y}{\delta} = 1} = {\beta \approx 1}} & \left( {{Equation}\quad 1} \right)\end{matrix}$The value of the constant, β, chosen will depend upon configuration andnumerical stability concerns. For the case of a flat plate, anappropriate value for this constant is 0.9975.

A turbulent boundary layer is characterized by a spectrum of pressureand shear fluctuations, the frequency, phase speed, and amplitudecharacteristics of which are a function of such factors as freestreamvelocity, body configuration, surface conditions, and pressure gradient.With flow over a rigid surface, there is no motion of the surface. Withan elastic or viscoelastic surface, the wall pressure and shearfluctuations act as a forcing function which can deform the surface,creating surface waves. With a viscoelastic surface, energy from theturbulent boundary layer may be absorbed and dissipated by the coating,thus necessitating proper specification of boundary conditions for bothReynolds stresses and energy absorption at the wall (i.e., kinematic anddynamic boundary conditions).

FIG. 1 is a schematic of a passive viscoelastic coating that interactswith a turbulent boundary layer of thickness, δ. For simplicity, thecase of a flow with a freestream velocity of U_(∞) over a coated flatplate is considered, where coordinates x, y, and z are in thelongitudinal, normal, and transverse directions, respectively.Alternatively, x is represented as x₁, y as x₂, and z as x₃. It isassumed that the fluid is viscous and incompressible, and that thecoating material is viscoelastic (i.e., possesses a combination ofelastic and viscous physical properties). If the surface is rigid, therewill be no oscillation or energy absorption. If the surface is elasticor viscoelastic, there will be surface oscillation and, if the surfaceis viscoelastic, there will be energy exchange as well. The componentsof longitudinal and vertical surface displacement of the viscoelasticmaterial are given by ξ₁ and ξ₂, and the friction velocity is definedas: $\begin{matrix}{u_{*} = \sqrt{\frac{\tau_{w}}{\rho}}} & \left( {{Equation}\quad 2} \right)\end{matrix}$where τ_(w) is the shear stress at the wall and ρ is the density of thefluid. The interaction of the turbulent flow with a viscoelastic coatingleads to the formation of a quasi-periodic surface wave. The motion andenergy absorption of the coating (kinematic and dynamic boundaryconditions) in turn affect the energy balance in the turbulent boundarylayer and the value of the friction drag, the latter of which is thesurface integral of the wall shear stress.

The methodology for this invention is schematically shown in FIG. 2. Thepresent invention includes solutions for: 1) turbulent boundary layer(TBL) parameters, including friction drag over rigid, viscoelastic, orelastic plates; and, 2) energy absorption and oscillation amplitudes ofa viscoelastic (VE) plate excited by a periodic load which approximatesthat of a turbulent boundary layer. These two parts of the solution arecoupled by boundary, conditions, parts of the solution are coupled byboundary conditions, both for energy absorbed by the surface and for theamplitudes of surface motion, and are part of an overall methodology forselecting drag-reducing coatings and for quantifying drag reduction forgiven flow conditions. Below, the two parts of the solution, as well asthe methodology that couples them, are described.

Characterization of a Turbulent Boundary Layer Over Rigid, Elastic, orViscoelastic Surfaces (TBL Problem)

In FIG. 2, TBL parameters are first characterized for a rigid surface(step 1) and, once material properties are defined, for a viscoelasticmaterial (step 3). The same basic approach is applied in step 3 as instep 1, except the kinematic and dynamic boundary conditions are defineddifferently. For a rigid surface, there is no surface motion and noenergy absorption. For an elastic surface, there is surface motion, butno energy absorption. For a viscoelastic surface, there is surfacemotion and energy absorption.

General System of Equations of Continuity, Motion and Energy: Turbulentflow parameters are obtained through the solution of a system ofequations of continuity, motion, and energy, with accompanying boundaryconditions. These equations are developed from principles ofconservation of mass, conservation of momentum (Newton's second law),and energy balance (as based on the first law of thermodynamics).

In Cartesian coordinates, the general equation of continuity for acompressible fluid having a density of ρ, and velocity components U, V,and W in the streamwise, normal, and transverse directions, is given byEquation 3, below. $\begin{matrix}{{\frac{\partial\rho}{\partial t} + \frac{\partial\left( {\rho\quad U} \right)}{\partial x} + \frac{\partial\left( {\rho\quad V} \right)}{\partial y} + \frac{\partial\left( {\rho\quad W} \right)}{\partial z}} = 0} & \left( {{Equation}\quad 3} \right)\end{matrix}$Alternatively, Equation 3 can be written in terms of indicial notation(as in Equation 4), where x, y, and z are represented by x₁, x₂, and x₃,respectively, an where U, V, and W are represented by U₁, U₂, and U₃,respectively. It is implied that the index, i, can have a value of 1, 2,or 3, and that a repeated index of i indicates summation.$\begin{matrix}{{\frac{\partial\rho}{\partial t} + \frac{\partial\left( {\rho\quad U_{i}} \right)}{\partial x_{i}}} = 0} & \left( {{Equation}\quad 4} \right)\end{matrix}$

Equation (3) is further simplified for a fluid that is incompressible,i.e., where the density of the fluid is constant, the following applies:$\begin{matrix}{\frac{\partial U_{i}}{\partial x_{i}} = 0} & \left( {{Equation}\quad 5} \right)\end{matrix}$

The generalized equations of motion, termed the Navier-Stokes equations,are expressed in Cartesian coordinates for the case of an incompressiblefluid with constant viscosity as: $\begin{matrix}{{\frac{\partial U_{i}}{\partial t} + {U_{j}\frac{\partial U_{i}}{\partial x_{j}}}} = {{{- \frac{1}{\rho}}\frac{\partial{??}}{\partial x_{i}}} + g_{i} + {v\frac{\partial^{2}U_{i}}{{\partial x_{j}}{\partial x_{j}}}}}} & \left( {{Equations}\quad 6a\text{-}6b\text{-}6c} \right)\end{matrix}$In Equations 6:

-   -   P is the mean pressure,    -   g_(i) is the body force vector, due to external fields, such as        gravity, which act on the element, and    -   ν is the kinematic viscosity (assumed to be constant).        Note that the above expression, written in indicial notation,        actually represents three equations, for the three components of        velocity in the x, y, and z directions.

Turbulent velocity components may be described as the sum of the meanand fluctuating components, U_(i) and u′, respectively, where U₁, U₂,and U₃, are equivalent to U, V, and W and where u′₁, u′₂, and u′₃ areequivalent to u′, v′, and w′:U_(i)=U_(i)+u_(i)  (Equations 7a-7b-7c)An overbar indicates time-averaging:U_(i)=Ū_(i)  (Equations 8a-8b-8c)Substituting Equations 7a-7c into Equations 6a-6c and time-averagingyields the following system of three complex nonlinear second-orderpartial differential equations of motion for turbulent flow:$\begin{matrix}{{U_{j}\frac{\partial U_{i}}{\partial x_{j}}} = {{{- \frac{1}{\rho}}\frac{\partial{??}}{\partial x_{i}}} + g_{i} + {v\frac{\partial^{2}U_{i}}{{\partial x_{j}}{\partial x_{j}}}} - \frac{\partial\overset{\_}{u_{i}^{\prime}u_{j}^{\prime}}}{\partial x_{j}}}} & \left( {{Equations}\quad 9a\text{-}9b\text{-}9c} \right)\end{matrix}$In Equations 9a-9c, the components u′² , v′² , and w′² are termed normalReynolds stresses, and the components in the form − u′v′, − v′w′, and −u′w′ are termed Reynolds shear stresses. The turbulent kinetic energy,k, is defined as: $\begin{matrix}{k = {\frac{1}{2}\left( {\overset{\_}{u^{\prime 2}} + \overset{\_}{v^{\prime 2}} + \overset{\_}{w^{\prime 2}}} \right)}} & \left( {{Equation}\quad 10} \right)\end{matrix}$Closure of the generalized system of equations including the continuityequation (Equation 5) and the equations of motion (Equations 9a-9c) fora turbulent flow requires seven additional equations to characterize thesix Reynolds stresses, u_(i)′u_(j)′, and the rate of transport of theturbulent kinetic energy, k. This invention solves for the isotropicdissipation rate, ε, which is related to energy transport through thefluid and at the fluid-surface interface. The energy transport equationis based upon the first law of thermodynamics, where heat, dQ, added toa volume during an element of time, dt, serves to increase internalenergy, dE, and to perform work, dWK. $\begin{matrix}{\frac{\mathbb{d}Q}{\underset{heat}{\underset{︸}{\mathbb{d}t}}} = {\frac{\mathbb{d}E}{\underset{energy}{\underset{︸}{\mathbb{d}t}}} + \frac{\mathbb{d}{WK}}{\underset{work}{\underset{︸}{\mathbb{d}t}}}}} & \left( {{Equation}\quad 11} \right)\end{matrix}$

There exist multiple approaches within the literature for developingadditional equations for Reynolds stress terms in turbulent flow, butthis invention adopts a Reynolds-stress-transport-type methodology. Inthis methodology, equations for Reynolds stresses take the followinggeneral form: $\begin{matrix}{{U_{k}\frac{\partial\overset{\_}{u_{i}^{\prime}u_{j}^{\prime}}}{\partial x_{k}}} = {P_{ij} - \Pi_{ij} - \frac{\partial J_{ijk}}{\partial x_{k}} - {2ɛ_{ij}}}} & \left( {{Equations}\quad 12a\text{-}12f} \right)\end{matrix}$where P_(ij) is termed the production, Π_(ij) is termed thepressure-strain correlation tensor, J_(ijk) is termed the diffusive fluxof the Reynolds stresses, and ε_(ij) is termed the dissipation tensor.

In the general case, equations for all six Reynolds stress terms, andfor the energy dissipation rate must expressed. The equation for theisotropic dissipation rate, ε, is similar in structure to the equationsfor the transport of Reynolds stresses. Full mathematical expressionsfor the Reynolds stress and isotropic dissipation rate equations shallbe expressed in the following section for the specific case of atwo-dimensional turbulent boundary layer.

In summary, the equations which are solved to determine turbulentboundary layer parameters include:

-   -   the continuity equation, Equation 5,    -   the three equations of motion, Equations 9a, 9b, 9c,    -   six equations for Reynolds shear and normal stresses, Equations        12a-12f, and    -   the equation for isotropic dissipation rate, ε (Equation 16        below).

The methodology for the solution of turbulent flow parameters involves afinite difference approximation of the system of equations of motion andcontinuity, with accompanying boundary conditions.

Turbulent Boundary Layer Equations: Complete mathematical formulationsare provided for the specific case of a turbulent boundary layer with asteady, two-dimensional mean flow and a constant freestream velocity,U_(∞). Two-dimensional turbulent boundary layer equations, as termed inthe literature, are derived from the general continuity equation(Equation 5) and equations of motion (Equations 9a-9c), given theassumptions that:

-   -   The mean transverse velocity, W, is zero.    -   Gravitational forces can be neglected.    -   The pressure gradient is approximately zero in the y direction.    -   The mean velocity in the streamwise direction, U, is much        greater than the mean velocity in the normal direction, V.    -   The rate of change of parameters in the x direction is much        smaller than the rate of change of parameters in the y        direction.        The above assumptions permit simplification of the set of        equations required to solve for turbulent boundary layer        parameters to include the revised continuity equation, Equation        (13): $\begin{matrix}        {{\frac{\partial U}{\partial x} + \frac{\partial V}{\partial y}} = 0} & \left( {{Equation}\quad 13} \right)        \end{matrix}$        and the equation of motion for the U velocity component        (Equation 14): $\begin{matrix}        {{{U\frac{\partial U}{\partial x}} + {V\frac{\partial U}{\partial y}}} = {{{- \frac{1}{\rho}}\frac{\partial{??}}{\partial x}} + {v\frac{\partial^{2}U}{\partial y^{2}}} - \frac{\partial\overset{\_}{u^{\prime}v^{\prime}}}{\partial y}}} & \left( {{Equation}\quad 14} \right)        \end{matrix}$

Where transport equations for the six Reynolds stress components arerequired in the general case, the Reynolds shear stress components −v′w′ and − u′w′ are considered to be very small, so that only equationsfor u′² , v′² , w′² , and − u′v′ are formulated in the format ofEquation 12 (repeated below as Equations 15a-15d). $\begin{matrix}{{U_{k}\frac{\partial\overset{\_}{u_{i}^{\prime}u_{j}^{\prime}}}{\partial x_{k}}} = {P_{ij} - \Pi_{ij} - \frac{\partial J_{ijk}}{\partial x_{k}} - {2ɛ_{ij}}}} & \left( {{Equations}\quad 15a\text{-}15d} \right)\end{matrix}$where P_(ij) is the production term, Π_(ij) is the pressure-straincorrelation tensor, J_(ijk) is the diffusive flux of the Reynoldsstresses, and ε_(ij) is the dissipation tensor. A fifth equation for εis: $\begin{matrix}{{{U\frac{\partial ɛ}{\partial x}} + {V\frac{\partial ɛ}{\partial y}}} = {{C_{s\quad 1}f_{1}\frac{ɛ}{k}P_{\Sigma}} - {C_{s\quad 2}f_{2}{\frac{ɛ}{k}\left\lbrack {ɛ - {v\frac{\partial^{2}k}{\partial y^{2}}}} \right\rbrack}} + {\frac{\partial}{\partial y}\left( {2ɛ_{t}\frac{\partial ɛ}{\partial y}} \right)} + {v\frac{\partial^{2}ɛ}{\partial y^{2}}}}} & \left( {{Equation}\quad 16} \right)\end{matrix}$where the expression for viscous diffusion may alternatively beapproximated as: $\begin{matrix}{{\left\lbrack \left\lbrack \quad \right.\quad \right.v\frac{\partial^{2}k}{\partial y^{2}}} = {2{v\left( \frac{\partial k^{1/2}}{\partial y} \right)}}} & {\left( {{Equation}\quad 17} \right)\left. \quad \right\rbrack\left. \quad \right\rbrack} \\{{v\frac{\partial^{2}k}{\partial y^{2}}} \approx {2{v\left( \frac{\partial k^{1/2}}{\partial y} \right)}}} & \left( {{Equation}\quad 17} \right)\end{matrix}$if required for numerical stability in solutions of viscoelastic,non-oscillating surfaces with limited grid points in the near-wallregion.

In Equations (15a-15d), the term P_(ij) may be expressed as:$\begin{matrix}{P_{ij} = {{{- \overset{\_}{u_{i}^{\prime}u_{k}^{\prime}}}\frac{\partial U_{j}}{\partial x_{k}}} - {\overset{\_}{u_{j}^{\prime}u_{k}^{\prime}}\frac{\partial U_{i}}{\partial x_{k}}}}} & \left( {{Equations}\quad 18a\text{-}18d} \right)\end{matrix}$In Equation (16), the term P_(Σ) may be expressed as: $\begin{matrix}{P_{\Sigma} = {\frac{1}{2}P_{ii}}} & \left( {{Equation}\quad 19} \right)\end{matrix}$

The pressure-strain correlation tensor, Π_(ij), which redistributesenergy between different components of Reynolds stresses, may beexpressed as: $\begin{matrix}{\Pi_{ij} = {{{C_{1}\left( \frac{ɛ}{k} \right)}\left( {\overset{\_}{u_{i}^{\prime}u_{j}^{\prime}} - {\delta_{ij}\frac{2}{3}k}} \right)} + {C_{2}\left( {P_{ij} - {\delta_{ij}\frac{2}{3}P_{\Sigma}}} \right)} + \pi_{{ij},1}^{\prime} + \pi_{{ij},2}^{\prime} + \pi_{{ij},3}^{\prime}}} & \left( {{Equations}\quad 20a\text{-}20d} \right)\end{matrix}$where the π′_(ij,1) terms represent near-wall redistribution ofturbulent energy from the streamwise component to the normal andtransverse components, the π′_(ij,2) terms represent near-wall variationof the Reynolds stress tensor component production, and the π′_(ij,3)terms represent near-wall redistribution of turbulent energyproportional to local vorticity: $\begin{matrix}\begin{matrix}{\pi_{{ij},1}^{\prime} = {{- C_{1}^{\prime}}\frac{ɛ}{k}\left( {{\overset{\_}{v^{\prime 2}}\delta_{ij}} - {\frac{3}{2}\left( {{\overset{\_}{v^{\prime}u_{i}^{\prime}}\delta_{jl}} + {\overset{\_}{v^{\prime}u_{j}^{\prime}}\delta_{il}}} \right)}} \right){f\left( \frac{l}{y} \right)}}} & \quad\end{matrix} & \left( {{Equations}\quad 21a\text{-}21d} \right) \\{\pi_{{ij},2}^{\prime} = {{- {C_{2}^{\prime}\left( {P_{ij} - {\frac{2}{3}\delta_{ij}P_{\Sigma}}} \right)}}{f\left( \frac{l}{y} \right)}}} & \left( {{Equations}\quad 22a\text{-}22d} \right) \\{{\left\lbrack \left\lbrack \quad \right.\quad \right.\pi_{{ij},3}^{\prime}} = {{- {C_{3}^{\prime}\left( {P_{ij} - D_{ij}} \right)}}f\quad\left( \frac{l}{y} \right)}} & {\left( {{Equations}\quad 23a\text{-}23d} \right)\left. \quad \right\rbrack\left. \quad \right\rbrack} \\{\pi_{{ij},3}^{\prime} = {{C_{3}^{\prime}\left( {P_{ij} - D_{ij}} \right)}{f\left( \frac{l}{y} \right)}}} & \left( {{Equations}\quad 23a\text{-}23d} \right)\end{matrix}$Here ${f\left( \frac{l}{y} \right)}\quad$is a unique damping function for the near-wall region: $\begin{matrix}{{f\left( \frac{l}{y} \right)} = {\frac{R_{t}}{R_{k}}\left( {1 + \sqrt{1 + \frac{100}{R_{t}}}} \right)}} & \left( {{Equation}\quad 24} \right)\end{matrix}$where: $\begin{matrix}{R_{k} = \frac{k^{1/2}y}{v}} & \left( {{Equation}\quad 25} \right)\end{matrix}$and $\begin{matrix}{R_{t} = \frac{k^{2}}{vɛ}} & \left( {{Equation}\quad 26} \right)\end{matrix}$Here, D_(ij) is a dissipation tensor: $\begin{matrix}\begin{matrix}{D_{ij} = {{{- \overset{\_}{u_{i}u_{j}}}\frac{\partial U_{i}}{\partial x_{j}}} - {\overset{\_}{u_{j}u_{i}}\frac{\partial U_{i}}{\partial x_{i}}}}} \\\frac{\partial J_{ijk}}{\partial x_{k}}\end{matrix} & \left( {{Equations}\quad 27a\text{-}27d} \right)\end{matrix}$is the gradient of turbulent and viscous diffusive flux of the Reynoldsstresses in the boundary layer, where only one component remains in theboundary-layer representation: $\begin{matrix}{{{- \frac{\partial J_{ijk}}{\partial x_{2}}} \pm \underset{{turbulent}\quad{diffusion}}{\underset{︸}{\frac{\partial}{\partial y}\left( {\left( {AC} \right)_{t}{\frac{k}{ɛ}\left\lbrack {\overset{\_}{v^{\prime 2}}\frac{\partial\overset{\_}{u_{i}^{\prime}u_{j}^{\prime}}}{\partial y}} \right\rbrack}} \right)}}} + \underset{{viscous}\quad{diffusion}}{v\underset{︸}{\frac{\partial\overset{\_}{u_{i}^{\prime}u_{j}^{\prime}}}{\partial y^{2}}}}} & \left( {{Equations}\quad 28a\text{-}28d} \right)\end{matrix}$where A is 6 in the equation for v′² , 2 for u′² and w′² , and 4 for −u′v′ (i.e., the effective gradients of turbulent diffusion are differentfor different components of Reynolds stress), and where the coefficientof turbulent diffusion, ε_(t), is: $\begin{matrix}{ɛ_{t} = {C_{t}\frac{k}{ɛ}\overset{\_}{v^{\prime 2}}}} & \left( {{Equation}\quad 29} \right)\end{matrix}$except for Equation (16), where:C_(i)=C_(ε)  (Equation 30)

The dissipation tensor, ε_(ij), is written as: $\begin{matrix}{ɛ_{ij} = {{v\overset{\_}{\frac{\partial\overset{\_}{u_{i}^{\prime}}}{\partial x_{k}}\frac{\partial\overset{\_}{u_{j}^{\prime}}}{\partial x_{k}}}} \approx {{f_{s}\frac{\overset{\_}{u_{i}^{\prime}u_{j}^{\prime}}}{2k}ɛ} + {\left( {1 - f_{s}} \right)\frac{1}{3}\delta_{ij}ɛ}}}} & \left( {{Equations}\quad 31a\text{-}31d} \right)\end{matrix}$where ƒ_(s) characterizes flow in the near-wall region: $\begin{matrix}{f_{s} = {\frac{1}{1 + {0.06R_{t}}}\quad{and}}} & \left( {{Equation}\quad 32} \right) \\{R_{t} = \frac{k^{2}}{vɛ}} & \left( {{Equation}\quad 33} \right)\end{matrix}$

Equation (16) includes two functions, ƒ₁ and ƒ₂, which also introducecorrections for near-wall flows:ƒ₁=1+0.8e^(−R) ^(t)   (Equation 34)ƒ₂=1−0.2e^(−R) ^(t) ²  (Equation 35)Values of constants for flow over a flat plate are as shown in Table 1:

C₁ C₂ C_(ε1) C_(a2) C_(t) C_(ε) C₁′ C₂′ C₃′ 1.34 0.8 1.45 1.9 0.12 0.150.36 0.45 0.036For the case of a two-dimensional boundary layer, turbulent boundarylayer parameters at different x and y locations are determined throughsolution of the continuity equation (Equation 13), the equation ofmotion in the x-direction (Equation 14), the transport equations for theReynolds stresses u′² , v′² , w′² , and − u′v′ (Equations 15a-15d), andthe equation for the isotropic energy dissipation rate (Equation 16),given appropriate boundary conditions. The problem is solved numericallyusing a finite difference approximation. All equations are reduced to astandard type of parabolic equation in terms of a given function, andsolution is obtained at designated grid points in an (x,y) coordinatesystem.

Boundary Conditions: Boundary conditions are values of parameters at thelimits of the boundary layer, i.e., at the surface and the freestream.The freestream velocity is defined as U_(∞). Boundary conditions at thesurface are specified for Reynolds normal and shear stresses (kinematicboundary conditions), as well as for the isotropic dissipation rate(dynamic boundary condition). For an arbitrary geometry, the x and ycoordinates of the surface must be specified. If the surface is a flatplate, the boundary will be along the line y=0.

Since oscillation amplitudes at the surface are small, linearizedkinematic boundary conditions, where mean velocities at the surface areassumed to be zero, are appropriate. Boundary conditions for fluctuatingvelocity components at the surface of a flat plate are expressed as:$\begin{matrix}{{u^{\prime}}_{y = 0} = {{\frac{\partial\xi_{1}}{\partial t}\cos\quad\Theta} - {\xi_{2}\frac{u_{*}^{2}}{v}}}} & \left( {{Equation}\quad 36} \right) \\{{v^{\prime}}_{y = 0} = \frac{\partial\xi_{2}}{\partial t}} & \left( {{Equation}\quad 37} \right) \\{{w^{\prime}}_{y = 0} = {\frac{\partial\xi_{1}}{\partial t}\sin\quad\Theta}} & \left( {{Equation}\quad 38} \right)\end{matrix}$where ξ₁ and ξ₂ are the longitudinal and vertical surface displacementcomponents, respectively, u* is the friction velocity (as previouslydefined), and Θ is the angle of the longitudinal axis relative to themean flow in the x₁-x₃ plane. With linearized boundary conditions, meanvelocities at the wall are assumed to be zero. Surface displacements areapproximated by the first mode of a Fourier series: $\begin{matrix}{\xi_{i} = {{\sum\limits_{j = 1}^{\infty}\quad{a_{ij}{\mathbb{e}}^{i\quad{\alpha_{j}{({x - {Ct}})}}}}} \approx {a_{il}{\mathbb{e}}^{i\quad{\alpha_{e}{({x - {Ct}})}}}}}} & \left( {{Equation}\quad 39} \right)\end{matrix}$Here, α₃ is the wavenumber corresponding to the maximum turbulent energyin the boundary layer, and is given by: $\begin{matrix}{\alpha_{e} = \frac{\omega_{e}}{C}} & \left( {{Equation}\quad 40} \right)\end{matrix}$where the energy-carrying frequency, ω_(e), is assumed to be:$\begin{matrix}{\omega_{e} = \frac{U_{\infty}}{\delta}} & \left( {{Equation}\quad 41} \right)\end{matrix}$and the phase speed corresponding to energy-carrying disturbances in theboundary layer is assumed to be:C≈0.8U_(∞)  (Equation 42)Since there is a range of frequencies which carry energy, as reportedwithin the scientific literature, it is advantageous to also performcalculations for the case where: $\begin{matrix}{\omega_{e} = {2\frac{U_{\infty}}{\delta}}} & \left( {{Equation}\quad 43} \right)\end{matrix}$In the absence of resonance, it is appropriate to time-averagecomponents of the Reynolds stress at the wall: $\begin{matrix}{{\overset{\_}{u^{\prime 2}}}_{y = 0} = {\frac{1}{2}{\omega_{e}^{2}\left( {{\xi_{1}}^{2} + {\frac{2u_{*}^{2}}{\omega_{e}v}{\xi_{1}}{\xi_{2}}{\sin\left( {\varphi_{2} - \varphi_{1}} \right)}} + {\frac{1}{\omega_{e}^{2}}\left( \frac{u_{*}^{2}}{v} \right)^{2}{\xi_{2}}^{2}}} \right)}}} & \left( {{Equation}\quad 44} \right) \\{{\overset{\_}{v^{\prime 2}}}_{y = 0} = {\frac{1}{2}\omega_{e}^{2}{\xi_{2}}^{2}}} & \left( {{Equation}\quad 45} \right) \\{{\overset{\_}{w^{\prime 2}}}_{y = 0} = {\frac{1}{2}\omega_{e}^{2}{\xi_{1}}^{2}\tan^{2}\Theta}} & \left( {{Equation}{\quad\quad}46} \right) \\{{{- \overset{\_}{u^{\prime}v^{\prime}}}}_{y = 0} = {\frac{1}{2}\omega_{e}^{2}{\xi_{2}}{\xi_{1}}{\cos\left( {\varphi_{2} - \varphi_{1}} \right)}}} & \left( {{Equation}\quad 47} \right)\end{matrix}$where |ξ_(i)| is the rms amplitude of the displacement. For a passiveisotropic viscoelastic coating excited by a forced load, the responsetakes the form of a traveling wave, so that the phase shift betweennormal and longitudinal displacements, φ₂-φ₁, will be approximately π/2,and the displacements, φ₂-φ₁ , will be approximately π/2, and theReynolds shear stresses at the surface will be approximately zero. Foranisotropic materials, the phase shift can be different, so thatnegative Reynolds shear stresses can be generated at the wall. For arigid wall, there will be no motion at the wall, so that Reynolds shearand normal stresses shall be equal to zero.

The boundary condition for the isotropic dissipation rate is:$\begin{matrix}{{{ɛ}_{y = 0} = {\frac{\partial}{\partial y}\left\lbrack {\left( {v + ɛ_{t}} \right)\frac{\partial k}{\partial y}} \right\rbrack}}}_{y = 0} & \left( {{Equation}\quad 48} \right)\end{matrix}$where the first term reflects viscous dissipation and the secondreflects absorption of energy by the viscoelastic material. For a rigidsurface, there is no energy absorption at the wall, so that the secondterm equals zero. The absorption of turbulent energy by the coating isequivalent to − p′v′, and can be approximated by${{{ɛ_{t}\frac{\partial k}{\partial y}}}_{\quad}}_{y = 0},$which is the diffusive flux of energy across the boundary, characterizedusing a gradient mechanism for turbulent diffusion. This expression ofthe dynamic boundary condition is compatible with the Reynolds stresstransport methodology of turbulence closure.

Equations 13, 14, 15a-15d, and 16 are solved for mean velocitycomponents, Reynolds normal and shear stresses, and energy dissipation,given the kinematic and dynamic boundary conditions (Equations (44)through (48)) based on the solution of the viscoelasticity problem (asdescribed in the following section). The problem is solved numerically,using finite difference approximations of the parabolic equations.Friction drag for a body with a viscoelastic coating is calculated asthe integral of wall shear stresses, τ_(w), over the surface of thebody, where: $\begin{matrix}{{{\tau_{w} = {\mu\frac{\partial U}{\partial y}}}}_{\quad}}_{y = 0} & \left( {{Equation}\quad 49} \right)\end{matrix}$for a two-dimensional body, and where μ=ρν is the dynamic viscosity.Comparison of results with those calculated for a rigid body ofidentical geometry under identical flow conditions leads to anestimation of anticipated friction drag reduction.

To reduce friction drag, it is necessary to minimize surface oscillationamplitudes, while maximizing the flux of turbulent energy from the flowinto the coating, − p′v′. If the amplitudes of surface oscillation, ξ₂⁺, do not exceed the thickness of the viscous sublayer, generally where:$\begin{matrix}{y^{+} = {\frac{{yu}_{*}}{v} < \approx 7}} & \left( {{Equation}\quad 50} \right)\end{matrix}$then the normal Reynolds stresses at the boundary (Equations (44)-(46))may be approximated as zero. For a coating which absorbs energy, withlow levels of oscillation, shear stresses in the near-wall region of theboundary layer decrease, as does the production of turbulence in theboundary layer. For a coating that oscillates at amplitudes greater thanthat of the viscous sublayer, the surface can act as a dynamic roughnesselement and thereby enhance the level of turbulence generated within theboundary layer.Response of a Viscoelastic Material to a Turbulent Boundary Layer (VEProblem)

The second part of the methodology determines the response of aviscoelastic material to a turbulent boundary layer (step 2 in FIG. 2).For a rigid surface, Reynolds stresses on the surface (Equations (44) to(47)) are zero, and the isotropic dissipation rate contains only theviscous term. However, for a viscoelastic material, the kinematic anddynamic boundary conditions are determined through solution of thetwo-dimensional conservation of momentum equation for a viscoelasticmaterial: $\begin{matrix}{{\rho_{s}\frac{\partial^{2}\xi_{i}}{\partial t^{2}}} = \frac{\partial\sigma_{ij}}{\partial x_{j}}} & \left( {{Equations}\quad 51a\text{-}51b} \right)\end{matrix}$where ρ_(s) is the material density, ξ₂ and ξ₂ are the longitudinalwhere ρ_(s) is the material density, ξ ₁ and ξ ₂ are the longitudinaland normal displacements through the thickness of the coating, andσ_(ij) is the amplitude of the stress tensor. The stress tensor for aviscoelastic material is written for a Kelvin-Voigt type of material as:σ_(ij)=λ(ω)ε^(s)δ_(ij)+2μ(ω)ε_(ij) ^(s)  (Equations 52a-52d)where ε_(ij) ^(s) is the strain tensor: $\begin{matrix}{ɛ_{ij}^{s} = {\frac{1}{2}\left( {\frac{\partial\xi_{i}}{\partial x_{j}} + \frac{\partial\xi_{j}}{\partial x_{i}}} \right)}} & \left( {{Equations}\quad 53a\text{-}53d} \right)\end{matrix}$and:ε^(s)=ε_(ii) ^(s)  (Equation 54)λ(ω) is the frequency-dependent Lame constant, which is defined in termsof the bulk modulus, K(ω), which can be reasonably approximated as thestatic bulk modulus, K₀, and the complex shear modulus, μ(ω):$\begin{matrix}{{\lambda(\omega)} = {K_{0} - {\frac{2}{3}{\mu(\omega)}}}} & \left( {{Equation}\quad 55} \right)\end{matrix}$

Displacements, ξ_(i), are approximated as periodic, in the form ofEquation (39), and can be expressed as a function of potentials oflongitudinal and transverse (shear) waves: $\begin{matrix}{\overset{\_}{\xi} = {{{\nabla\phi} + {{\nabla\overset{\rightarrow}{\Psi}}\overset{\rightarrow}{\Psi}}} = \left\{ {0,0,\Psi} \right\}}} & \left( {{Equation}\quad 56} \right)\end{matrix}$where ∇_(φ) is the gradient of φ and ∇×{right arrow over (Ψ)} is thecurl of the vector {right arrow over (Ψ)}. Equation (54) can berewritten as two decoupled equations for the two wave potentials:$\begin{matrix}{{\left\lbrack {{2{\mu(\omega)}} + {\lambda(\omega)}} \right\rbrack{\nabla^{2}\varphi}} = {\rho_{s}\frac{\partial^{2}\varphi}{\partial t^{2}}}} & \left( {{Equation}\quad 57} \right) \\{{{\mu(\omega)}{\nabla^{2}\Psi}} = {\rho_{s}\frac{\partial^{2}\Psi}{\partial t^{2}}}} & \left( {{Equation}\quad 58} \right)\end{matrix}$

Equations (57) and (58) can be solved for the potentials, φ and Ψ, andhence for displacements, velocities, and stresses through the thicknessof the coating, if boundary conditions are specified. The coating isfixed at its base, so that the longitudinal and normal displacements arezero, and the shear stress and pressure load on the surface is known.Pressure and shear pulsations on the coatings are approximated asperiodic functions, with a form similar to that of the displacements inEquation (39), but with the following magnitudes, respectively:τ_(ω)=ρu*²  (Equation 59)ρ_(rms)=K_(p)τ_(ω-)K_(p)ρu*²  (Equation 60)I_(ω) =ρμ _(φ) ²  (Equation 59 )ρ_(rms) =K _(p) I _((t)) =K _(p) ρμ* ²  (Equation 60 )where K_(p) is the Kraichnan parameter, whose value is approximated as2.5.

If shear pulsations are included, a phase shift between shear andpressure pulsations must also be introduced.

If calculations are performed for a unit load, then surfacedisplacements under actual load will be: $\begin{matrix}{{{\xi_{i}}❘_{p_{{rms} = {actual}}}} = {{\frac{{\rho\quad K_{p}u_{*}^{2}{\xi_{i}}}❘_{p_{rms} = 1}}{H}H{\overset{\sim}{u}}_{*}^{2}} = {C_{ki}H{\overset{\sim}{u}}_{*}^{2}}}} & \left( {{Equation}\quad 61} \right)\end{matrix}$Kinematic boundary conditions in Equations (44) to (47) for theturbulent boundary layer problem are rewritten in terms of output fromthe materials problem: $\begin{matrix}{\frac{\overset{\_}{u^{\prime 2}}}{U_{\infty}^{2}} = {\frac{1}{2}\left( \frac{H}{\delta} \right)^{2}\left( {C_{k1} + {{\overset{\sim}{u}}_{*}{Re}_{*}C_{k2}}} \right)^{2}{{\overset{\sim}{u}}_{*}^{4}\left( {{Re}_{*} = \frac{u_{*}\delta}{v}} \right)}}} & \left( {{Equation}\quad 62} \right) \\{\frac{\overset{\_}{v^{\prime 2}}}{U_{\infty}^{2}} = {\frac{1}{2}\left( \frac{H}{\delta} \right)^{2}C_{k2}^{2}{\overset{\sim}{u}}_{*}^{4}}} & \left( {{Equation}\quad 63} \right) \\{\frac{\overset{\_}{w^{\prime 2}}}{U_{\infty}^{2}} = {\frac{1}{2}\left( \frac{H}{\delta} \right)^{2}C_{k1}^{2}{\overset{\sim}{u}}_{*}^{4}\tan^{2}\Theta}} & \left( {{Equation}\quad 64} \right) \\{{- \frac{\overset{\_}{u^{\prime}v^{\prime}}}{U_{\infty}^{2}}} = 0} & \left( {{Equation}\quad 65} \right)\end{matrix}$

Dynamic boundary conditions are rewritten in the form: $\begin{matrix}{{{{- \frac{{\overset{\_}{p^{\prime}v^{\prime}}}_{y = 0}}{{\rho U}_{\infty}^{3}}} \approx {{\overset{\sim}{ɛ}}_{t}\frac{\partial\overset{\sim}{k}}{\partial\overset{\sim}{y}}}}}_{\quad}}_{y = 0} & \left( {{Equation}\quad 66} \right)\end{matrix}$where:C_(k3)=C_(k2)K_(p)γ(ω)and where γ(ω) is a dissipative function of the coating material.

The flux of turbulent fluctuating energy through the surface can besolved for directly, as:− p′v′|_(y=0)=−¼|p′|ω(−tξ₂+tξ*₂)  (Equation 67)but the nondimensionalized flux can also be approximated as a diffusiveflux term, using the gradient diffusion approach: $\begin{matrix}{{{{- \frac{{\overset{\_}{p^{\prime}v^{\prime}}}_{y = 0}}{{\rho U}_{\infty}^{3}}} \approx {{\overset{\sim}{ɛ}}_{t}\frac{\partial\overset{\sim}{k}}{\partial\overset{\sim}{y}}}}}_{\quad}}_{y = 0} & \left( {{Equation}\quad 68} \right)\end{matrix}$where: $\begin{matrix}{\overset{\sim}{k} = \frac{k}{U_{\infty}^{2}}} & \left( {{Equation}\quad 69} \right) \\{\overset{\sim}{y} = \frac{y}{\delta}} & \left( {{Equation}\quad 70} \right) \\{{\overset{\sim}{ɛ}}_{t} = \frac{ɛ_{t}}{U_{\infty}\delta}} & \left( {{Equation}\quad 71} \right)\end{matrix}$

Equation 68 provides a basis for determining the value of the kinematiccoefficient of turbulence diffusion, ε_(q), on an absorbing surface.Substituting Equation 68 into Equation 66 yields the followingexpression for ε_(q), defined as ε_(t)|_(y=0): 66 yields the followingexpression for {tilde over (ε)}_(q) , defined as {tilde over (ε)}_(t)|_(y=0): $\begin{matrix}{{\overset{\sim}{ɛ}}_{q} = {\frac{ɛ_{t}❘_{y = 0}}{U_{\infty}\delta} = {\frac{1}{4}{C_{k3}\left( \frac{H}{\delta} \right)}^{2}{\overset{\sim}{u}}_{*}\frac{{\overset{\sim}{u}}_{*}^{3}\Delta_{\max}}{k_{\max}^{+} - k_{q}^{+}}}}} & \left( {{Equation}\quad 72} \right)\end{matrix}$Here, K⁺ _(max) is the maximum of turbulence kinetic energy, and K⁺ _(q)is the kinetic energy of the oscillating surface, both quantitiesnondimensionalized by U_(∞) ². y_(max) ⁺ is defined as the normaldistance from the surface to the maximum of turbulence energy,nondimensionalized as follows: $\begin{matrix}{y_{\max}^{+} = {\frac{y_{\max}}{v}u_{*}}} & \left( {{Equation}\quad 73} \right)\end{matrix}$and: $\begin{matrix}{\Delta_{\max} = \frac{y^{+}}{{Re}_{*}}} & \left( {{Equation}\quad 74} \right)\end{matrix}$where: $\begin{matrix}{{Re}_{*} = \frac{u_{*}\delta}{v}} & \left( {{Equation}\quad 75} \right)\end{matrix}$Thus we can determine the dissipation rate at the wall based on(Equation 48).Methodology to Choose Properties of a Drag-Reducing ViscoelasticMaterial

A methodology to choose properties of a viscoelastic coating thatreduces turbulent friction drag necessarily requires both of thepreviously described solutions for turbulent boundary layer parametersand response of a viscoelastic material.

For the case of two-dimensional flow over a flat plate, the TBL problemis solved for a rigid plate, in order to determine the boundary layerthickness, δ, at a given freestream velocity, U_(∞), and location. Theboundary layer thickness is determined from the finite differencesolution of the seven equations of continuity, motion in thex-direction, transport equations for Reynolds normal and shear stresses,and the equation for dissipation rate, assuming no motion at the wall.The external limit of the boundary layer is defined as that locationwhere the ratio of the mean velocity to the freestream velocity is aconstant, β, between 0.95 and 1.0.

The frequency-dependent, complex shear modulus of a material, μ(ω), canbe expressed in different mathematical forms, some of which approximateexperimentally measured shear modulus data more accurately than others.A single relaxation time (SRT) material is one where the complex shearmodulus is expressed using a single relaxation time, τ_(s), and a singlevalue for the dynamic shear modulus, μ₂. In and a single value for thedynamic shear modulus μ₂ . In Equation (76), an SRT material would berepresented for the case of N=1. A multiple relaxation time (MRT)material is one where N>1 in the representation for complex shearmodulus in Equation (78). $\begin{matrix}{{\mu(\omega)} = {\mu_{0} + {\sum\limits_{j = 1}^{N}\quad{\mu_{j}\left\lbrack {\frac{\left( {\omega\tau}_{j} \right)^{2}}{1 + \left( {\omega\tau}_{j} \right)^{2}} + {{\mathbb{i}}\frac{{\omega\tau}_{j}}{1 + \left( {\omega\tau}_{j} \right)^{2}}}} \right\rbrack}}}} & \left( {{Equation}\quad 76} \right)\end{matrix}$The Havriliak-Negami (HN) representation for the complex shear modulus,is given by Equation 74. This equation is more complex, but often ismore suitable for describing real materials: $\begin{matrix}{\frac{\mu - \mu_{\infty}}{\mu_{0} - \mu_{\infty}} = \frac{1}{\left\lbrack {1 + \left( {\mathbb{i}\omega\tau}_{HN} \right)^{\alpha_{HN}}} \right\rbrack^{\beta_{HN}}}} & \left( {{Equation}\quad 77} \right)\end{matrix}$For a HN type of material, whose complex shear modulus is expressed inthe form of Equation (77), μ_(∞) is the limiting high-frequency modulusand α_(HN) and β_(HN) are constants. For the type of polymeric materialsused for drag-reducing coatings, K(ω) is essentially constant, with avalue of approximately 1×10⁸ Pa.

It is recommended to first determine an optimal SRT type of material,and then to choose an HN type of material, whose properties can becreated with available polymer chemistry.

An SRT material can be adequately characterized by the materialthickness, H, the density, ρ_(s), the static shear modulus, μ₀, thedynamic shear modulus, μ_(s), and the relaxation time, τ. An appropriatedensity for the viscoelastic material, ρ_(s), is within 10% that ofwater. For an SRT material, the initial guess for the static shearmodulus of the material, μ₀, is:μ₀=ρ_(s) C ²  (Equation 78)based on the criterion that the speed of shear waves in the material isapproximately the same as the phase speed of the energy-carryingdisturbances, C. This phase speed, C, is assumed to be 0.8 of the valueof the freestream velocity, U_(∞) (Equation 42). If the convectivevelocity exceeds the shear wave velocity, an instability occurs, andlarge waves appear on the surface of the material, leading to anincrease of drag for the coating.

An initial choice for thickness, H, for isotropic viscous materials,where $\frac{\mu_{s}}{\mu_{0}} > 1$is: $\begin{matrix}{H = \frac{3C}{\omega_{e}}} & \left( {{Equation}\quad 79} \right)\end{matrix}$and for isotropic, low viscosity materials, where$\frac{\mu_{s}}{\mu_{0}} < 1$is: $\begin{matrix}{H = \frac{5C}{\omega_{e}}} & \left( {{Equation}\quad 80} \right)\end{matrix}$

The optimal desired thickness for a coating may be greater thanpractical for a given application. While isotropic coatings thinner thanrecommended in Equations 79-80 can still be effective, anisotropiccoatings that are stiffer in the normal dimension relative to thetransverse and longitudinal dimensions can provide equivalentperformance with significant reduction in thickness.

Given specified values of H, μ₀, and ρ_(s), the VE problem, as expressedin Equations (57) and (58), is solved numerically for a matrix of valuesof τ_(s) and μ_(s) (i.e., for different values of the complex shearmodulus) and for a range of wavenumbers. The wavenumber corresponding tothe maximum turbulent energy in the boundary layer is: $\begin{matrix}{\alpha_{e} = \frac{\omega_{e}}{C}} & \left( {{Equation}\quad 81} \right)\end{matrix}$where the frequency, ω_(e), for maximum energy-carrying disturbances isestimated by Equation (41). Calculations yield surface displacementamplitudes and the flux of turbulent fluctuating energy into thecoating. The best combination of properties for a SRT material occurswhere the surface displacement under actual load (Equation (61)) is lessthan the viscous sublayer thickness, and where the energy flux into thecoating (Equation (57)) is at a maximum. Furthermore, it is desirable tomaintain this criterion for a range of frequencies from approximatelyone decade below to one decade above the energy-carrying frequency,ω_(e).

Once a set of optimal values of τ_(s) and μ_(s) are determined for a setof specified values of H, μ₀, and ρ_(s), the calculations are iteratedusing slightly different values of thickness, H, and static modulus, μ₀.From these calculations are chosen the optimal set of parameters for anSRT material (H, μ₀, ρ_(s), τ_(s) and μ_(s)), given specified flowconditions and configuration.

The complex shear moduli of real polymeric materials, such aspolyurethanes and silicones, which are candidates for viscoelasticcoatings cannot be adequately described by the SRT representation. Morecomplex MRT or HN representations of the shear modulus require multipleconstants, and are less suitable for numerical parametric evaluation.Therefore, results for SRT materials are used to select candidatematerials, such as described by the HN formulation (Equation (79)) whichcan be more readily fabricated in practice. As a guideline, it isdesired to match the complex shear modulus curves of the target SRTmaterial and the HN material (value and slope) over frequencies rangingfrom decade below to one decade above ω_(e), with the most importantmatching being in the immediate vicinity of ω_(e).

To design a multi-layer isotropic coating, properties of a complex shearmodulus, density, and thickness are specified for individual layers, andnon-slip boundary conditions between layers are imposed. The propertiesof the upper layer are specified according to the methodology for asingle layer, and the lower layers will have progressively lower staticshear moduli, as optimized for lower freestream velocities. Thus,well-designed multi-layer coatings can reduce drag over a range offreestream velocities.

In the design of an anisotropic coating, the complex shear modulus hasdifferent values in the normal direction relative to the longitudinaland transverse directions (hereinafter termed transversely isotropic).If the viscoelastic material follows a single-relaxation time model,then the static shear modulus, μ₀, the dynamic shear modulus, μ_(s), andthe relaxation time, τ_(s), will differ with direction, as expressed inEquations (82) and (83). The static shear modulus in the normaldirection, μ₀₁, will be greater than than in the longitudinal-transverseplane, μ₀₂. The complex shear modulus in the normal direction isexpressed as: $\begin{matrix}{{\mu_{1}(\omega)} = {\mu_{01} + {\mu_{s1}\left\lbrack {\frac{\left( {\omega\tau}_{s1} \right)^{2}}{1 + \left( {\omega\tau}_{s1} \right)^{2}} + {{\mathbb{i}}\frac{{\omega\tau}_{s1}}{1 + \left( {\omega\tau}_{s1} \right)^{2}}}} \right\rbrack}}} & \left( {{Equation}\quad 82} \right)\end{matrix}$while the shear modulus in the streamwise and transverse directions isexpressed as: $\begin{matrix}{{\mu_{2}(\omega)} = {\mu_{02} + {\mu_{s2}\left\lbrack {\frac{\left( {\omega\tau}_{s2} \right)^{2}}{1 + \left( {\omega\tau}_{s2} \right)^{2}} + {{\mathbb{i}}\frac{{\omega\tau}_{s2}}{1 + \left( {\omega\tau}_{s2} \right)^{2}}}} \right\rbrack}}} & \left( {{Equation}\quad 83} \right)\end{matrix}$For a viscoelastic, transversely isotropic material, surface oscillationamplitudes can be reduced relative to an isotropic material, while thelevel of energy flux into the material is increased. Thus, well-designedanisotropic coatings will be significantly thinner than isotropiccoatings associated with the same level of drag reduction.Methodology to Choose Structure of a Drag-Reducing Viscoelastic Coating

A further aspect of coating design is the choice of internal structurewithin the viscoelastic material. In practical applications ofviscoelastic coatings, the coating will be finite in length, withleading, trailing, and side edges. The influence of the finite edgesaffects coating performance. Well posed edges can order and stabilizetransverse and longitudinal vortical structures in the near-wall regionof the flow and thereby delay the deformation of these vorticalstructures and enhance the stability of the flow. However, unstructurededges can accentuate the amplitude of oscillations of the viscoelasticmaterial in this region. Local instabilities can degrade the performanceof the coating, so that, even with a well-designed material, theinfluence of the edges can lead to a drag increase. Hence, the coatingis structured in the vicinity of finite edges. The thickness of thecoating is decreased to minimize such oscillations, using techniquessuch as a rigid wedge underneath the coating, or other localizedstructure near an edge (FIG. 3). For large bodies, a continuous coatingmay be impractical or difficult to fabricate. An alternative design is apiecewise continuous coating, composed of finite segments of coating,where both the longitudinal and transverse edges of the coating systemare organized to stabilize flow structures and to minimize adverseeffects at the edges of each segment.

In addition to well-posed edges, viscoelastic coatings may be combinedwith surface structure to enhance the stabilization of longitudinalvortices along the length of the coating, and hence to increase thelevel of drag reduction through multiple physical mechanisms. Structurecan include the placement of riblets on top of the viscoelastic coating,or the creation of so-called “inverse” riblets. In the latter case, aviscoelastic coating may be molded over ribs or ridges of rigidmaterial, so that longitudinal riblet structures form when fluid flowsover the viscoelastic surface.

The dimensions (scales) of the segments and the dimensions of thestructures within the coating are selected as multiples of thetransverse and longitudinal scales in the near wall turbulent flow.These scales vary with body speed, position along the body and whennon-Newtonian additives, such as dilute aqueous solutions ofhigh-molecular weight polymers, are present.

1. A method for estimating the reduction in friction drag of a body witha viscoelastic coating as compared to the friction drag of a rigidsurface of identical size and shape as said body, during turbulent flowat a specified freestream velocity, said method comprising the followingsteps, performed in the order indicated: (I) using boundary conditionsfor a rigid surface, determining characteristics of a turbulent boundarylayer at the specified freestream velocity, said characteristicsincluding the boundary layer thickness, phase speed and frequencycorresponding to maximum-energy-carrying disturbances, mean velocityprofiles, Reynolds stress distributions, wall shear stress distribution,and friction drag, (II) selecting a density, complex shear modulus, andthickness of a viscoelastic coating which corresponds to minimumoscillation amplitudes and maximum flux of energy into the viscoelasticcoating during excitation by a forcing function substantially identicalto the excitation produced by the turbulent boundary layer as determinedin step I, and (III) using oscillation amplitudes and energy fluxcorresponding to the viscoelastic coating of selected density, complexshear modulus and thickness as determined in step II, determiningcharacteristics of the turbulent boundary layer at the specifiedfreestream velocity, including mean velocity profiles, Reynolds stressdistributions, wall shear stress distribution, and friction drag, and(IV) determining a percentage reduction in friction drag as the ratio ofi) the friction drag as determined in step I minus the friction drag asdetermined in step III, divided by ii) the friction drag as determinedin step I.
 2. The method of claim 1, wherein step I thereof includes thefollowing substeps, not necessarily performed in the order indicated:(a) for a rigid surface of specified geometry, and for a givenfreestream velocity, U_(∞), solving the equation of continuity:$\frac{\partial U_{i}}{\partial x_{i}} = 0$  and the equations of motionfor an incompressible fluid and steady flow with constant kinematicviscosity, ν, constant density, ρ, negligible body forces, and gradientsof pressure, P:${U_{j}\frac{\partial U_{i}}{\partial x_{j}}} = {{{- \frac{1}{\rho}}\frac{\partial{??}}{\partial x_{i}}} + g_{i} + {v\frac{\partial^{2}U_{i}}{{\partial x_{j}}{\partial x_{j}}}}}$where U_(i) refers to velocity components in the x, y, and z directions(x₁, x₂, and x₃ in indicial notation) and where fluid velocities at thesurface of the body are zero, (b) from the velocity field, U_(i),determined from the solution of the general continuity and motionequations in step 2(a), determining the boundary layer thickness as afunction of body geometry, where the edge of the boundary layer isdefined as that location where the ratio of mean velocity to thefreestream velocity is a constant, β, between 0.95 and 1.0:${\frac{U}{U_{\infty}}}_{\frac{y}{\delta} = 1} = \beta$ (c) from thesolution of the general continuity and motion equations in step 2(a),determining the shear stress, τ_(w), along the body, the localcoefficient of friction drag, and the integral drag, and (d) estimatingthe frequency, ω_(e), and the phase speed, C, corresponding tomaximum-energy-carrying disturbances in the boundary layer.
 3. Themethod of claim 2, wherein in substep 2(b), the boundary layer thicknessis approximated as that location where the ratio of the mean velocity tothe freestream velocity has a constant value of β=9975. velocity has aconstant value of β=0.9975.
 4. The method of claim 2, wherein in substep2(d), the maximum energy-carrying frequency for disturbances in theboundary layer is estimated as: $\omega_{e} = \frac{U_{\infty}}{\delta}$and the phase speed corresponding to energy-carrying disturbances in theboundary layer is assumed to be:C=0.8U_(∞).
 5. The method of claim (2), for the case of two-dimensionalflow over a flat plate, wherein the equations of continuity and theequations of motion are reduced to the following system of equations,where a Reynolds transport approach is adopted for turbulence closure,and where P is the mean pressure, ν is the kinematic viscosity, ρ is thedensity, U is the mean longitudinal velocity component, V is the meannormal velocity component, and u′, v′, and w′ are components offluctuating velocity in the streamwise, normal, and transversedirections, respectively, and ε is the isotropic dissipation rate:$\begin{matrix}{{\frac{\partial U}{\partial x} + \frac{\partial V}{\partial y}} = 0} & \quad \\{{{U\frac{\partial U}{\partial x}} + {V\frac{\partial U}{\partial y}}} = {{{- \frac{1}{\rho}}\frac{\partial{??}}{\partial x_{i}}} + {v\frac{\partial^{2}U}{\partial y^{2}}} - \frac{\partial\overset{\_}{u^{\prime}v^{\prime}}}{\partial y}}} & \quad \\{{U_{k}\frac{\partial\overset{\_}{u_{i}^{\prime}u_{j}^{\prime}}}{\partial x_{k}}} = {P_{ij} - \Pi_{ij} - \frac{\partial J_{ijk}}{\partial x_{k}} - {2ɛ_{ij}}}} & \quad \\{{{U\frac{\partial ɛ}{\partial x}} + {V\frac{\partial ɛ}{\partial y}}} = {{C_{ɛ\quad 1}f_{1}\frac{ɛ}{k}P_{\Sigma}} - {C_{ɛ\quad 2}f_{2}{\frac{ɛ}{k}\left\lbrack {ɛ - {v\frac{\partial^{2}k}{\partial y^{2}}}} \right\rbrack}} + {\frac{\partial}{\partial y}\left( {2ɛ_{t}\frac{\partial ɛ}{\partial y}} \right)} + {v\frac{\partial^{2}ɛ}{\partial y^{2}}}}} & \quad\end{matrix}$ where: P_(ij) is the production tensor, expressed as${\left\lbrack \left\lbrack \quad \right.\quad \right.P_{ij}} = {{{{- \overset{\_}{u_{i}^{\prime}u_{k}^{\prime}}}\frac{\partial U_{j}}{\partial x_{k}}} - {\overset{\_}{u_{j}^{\prime}u_{k}^{\prime}}\frac{\partial U_{i}}{\partial x_{k}}P_{\Sigma}}} = {\frac{1}{2}P_{ii}\left. \quad\left. \quad \right\rbrack \right\rbrack}}$$P_{if} = {{{- \overset{\_}{u_{i}^{\prime}u_{k}^{\prime}}}\frac{\partial U_{j}}{\partial x_{k}}} - {\overset{\_}{u_{j}^{\prime}u_{k}^{\prime}}\frac{\partial U_{i}}{\partial x_{k}}}}$$P_{\Sigma} = {\frac{1}{2}P_{ii}}$ Π_(ij) is the pressure-straincorrelation tensor, expressed as:$\Pi_{ij} = {{{C_{1}\left( \frac{ɛ}{k} \right)}\left( {\overset{\_}{u_{i}^{\prime}u_{j}^{\prime}} - {\delta_{ij}\frac{2}{3}k}} \right)} + {C_{2}\left( {P_{ij} - {\delta_{ij}\frac{2}{3}P_{\Sigma}}} \right)} - {\left( {{C_{i}^{\prime}\frac{ɛ}{k}\left( {{\overset{\_}{v^{\prime 2}}\delta_{ij}} - {\frac{3}{2}\left( {\overset{\_}{v^{\prime}u_{i}^{\prime}}\delta_{ij}} \right)}} \right)} + {C_{2}^{\prime}\left( {P_{ij} - {\frac{2}{3}\delta_{ij}P_{\Sigma}}} \right)} - {C_{2}^{\prime}\left( {P_{ij} - D_{ij}} \right)}} \right)\frac{R_{t}}{R_{k}}\left( {1 + \sqrt{1 + \frac{100}{R_{t}}}} \right)}}$where: k is the turbulent kinetic energy, expressed as$k = {\frac{1}{2}\left( {\overset{\_}{u^{\prime 2}} + \overset{\_}{v^{\prime 2}} + \overset{\_}{w^{\prime 2}}} \right)}$D_(ij) is the dissipation tensor, expressed as $\begin{matrix}{D_{ij} = {{{- \overset{\_}{u_{i}u_{j}}}\frac{\partial U_{i}}{\partial x_{j}}} - {\overset{\_}{u_{j}u_{i}}\frac{\partial U_{i}}{\partial x_{i}}}}} \\{R_{t} = \frac{k^{2}}{vɛ}}\end{matrix}$  is a nondimensional ratio between the square of thekinetic energy and the product of the kinematic viscosity anddissipation velocity $R_{k} = \frac{k^{1/2}y}{v}$  is the distance fromthe wall, nondimensionalized by the ratio of the kinematic viscosity tothe square of the kinetic energy C₁, C₂, C₁′, C₂′, and C₃′ areconstants, defined as C₁=1.34, C₂=0.8, C₁′=0.36, C₂′=0.45, and C₃′=0.036J_(ijk) is the tensor that describes the diffusive flux of the Reynoldsstresses, components of which are expressed as:${- \frac{\partial J_{ijk}}{\partial x_{1}}} \approx {0 - \frac{\partial J_{ijk}}{\partial x_{3}}} \approx {0 - \frac{\partial J_{ijk}}{\partial x_{2}}} \approx {\underset{{turbulent}\quad{diffusion}}{\underset{︸}{\frac{\partial}{\partial y}\left( {\left( {AC} \right)_{t}{\frac{\kappa}{ɛ}\left\lbrack {\overset{\_}{v^{\prime 2}}\frac{\partial\overset{\_}{u_{i}^{\prime}u_{j}^{\prime}}}{\partial y}} \right\rbrack}} \right)}} + \underset{{viscous}\quad{diffusion}}{v\underset{︸}{\frac{\partial^{2}\overset{\_}{u_{i}^{\prime}u_{j}^{\prime}}}{\partial y^{2}}}}}$where A is 6 for v′² , 2 for u′² and w′² , and 4 for − u′v′, where C₁ isa constant, defined as 0.12 ε_(ij) is the dissipation tensor, expressedas:$ɛ_{ij} = {{\frac{1}{1 + {0.06R_{t}}}\overset{\_}{\frac{u_{i}^{\prime}u_{j}^{\prime}}{2k}}ɛ} + {\left( {1 - \frac{1}{1 + {0.06R_{t}}}} \right)\frac{1}{3}\delta_{ij}ɛ}}$C_(ε1) and C_(ε2) are constants, defined as C_(ε1)=1.45 and C_(ε2)=1.9ε_(t) is the coefficient of turbulent diffusion, defined as${ɛ_{t} = {C_{t}\frac{k}{ɛ}\overset{\_}{v^{\prime 2}}}},$  except in theequation for the isotropic dissipation rate, where${ɛ_{t} = {C_{ɛ}\frac{k}{ɛ}\overset{\_}{v^{\prime 2}}}},$  and whereC_(ε)=0.15, ƒ₁=1+0.8e^(−R) ^(t) ƒ₂=1−0.2e^(−R) ^(t) ².
 6. The method ofclaim 5 wherein Reynolds stress boundary conditions are derived fromvalues of surface amplitudes, ξ₁ and ξ₂, calculated from solution of themomentum equation for a viscoelastic material: $\begin{matrix}{{\overset{\_}{u^{\prime 2}}}_{y = 0} = {\frac{1}{2}{\omega_{e}^{2}\left( {{\xi_{1}}^{2} + {\frac{2u_{*}^{2}}{\omega_{e}v}{\xi_{1}}{\xi_{2}}{\sin\left( {\varphi_{2} - \varphi_{1}} \right)}} + {\frac{1}{\omega_{e}^{2}}\left( \frac{u_{*}^{2}}{v} \right)^{2}{\xi_{2}}^{2}}} \right)}}} \\{{\overset{\_}{v^{\prime 2}}}_{y = 0} = {\frac{1}{2}\omega_{e}^{2}{\xi_{2}}^{2}}} \\{{\overset{\_}{w^{\prime 2}}}_{y = 0} = {\frac{1}{2}\omega_{e}^{2}{\xi_{1}}^{2}\tan^{2}\Theta}} \\{{{- \overset{\_}{u^{\prime}v^{\prime}}}}_{y = 0} = {\frac{1}{2}\omega_{e}^{2}{\xi_{2}}{\xi_{1}}{\cos\left( {\varphi_{2} - \varphi_{1}} \right)}}}\end{matrix}$ where: $\omega_{e} = \frac{U_{\infty}}{\delta}$ (φ₂-φ₁, isthe phase difference between normal and longitudinal oscillationamplitudes Θ is the angle between the mean flow and the x-axis u* is thefriction velocity, defined as μ_(φ) is the friction velocity, defined as${u_{*} = \sqrt{\frac{\tau_{w}}{\rho}}},$ where: ρ is the density of thefluid, and τ_(w) is the shear stress at the wall, expressed as${{\tau_{w} = {\rho\quad v\frac{\partial U}{\partial y}}}}_{y = 0}$|ξ_(i)| is the root-mean-square (rms) amplitude of the displacement. 7.The method of claim 6 wherein the energy boundary condition is expressedin terms of the solution of the momentum equation for a viscoelasticmaterial where: (a) the energy flux into the coating is determined fromthe pressure velocity correlation, which, given a periodic forcingfunction can be expressed in terms of the amplitude of the normaloscillation amplitude (and its complex conjugate, as designated by *):${{- \overset{\_}{p^{\prime}v^{\prime}}}}_{y = 0} = {{- \frac{1}{4}}{p^{\prime}}{\omega\left( {{- {\mathbb{i}\xi}_{2}} + {\mathbb{i}\xi}_{2}^{*}} \right)}}$(b) the amplitude of surface oscillations, ξ₂, given surface loadingcorresponding to the turbulent boundary layer, is expressed as:${{\xi_{2}}❘_{p_{{rms} = {actual}}}} = {{\frac{{\rho\quad K_{p}u_{*}^{2}{\xi_{2}}}❘_{p_{rms} = 1}}{H}H{\overset{\sim}{u}}_{*}^{2}} = {C_{k2}H{\overset{\sim}{u}}_{*}^{2}}}$ where: H is the thickness of the coating,${{\overset{\sim}{u}}_{*} = \frac{u_{*}}{U_{\infty}}},{and}$ u* is thefriction velocity, defined asμ₀ is the friction velocity, defined as${u_{*} = \sqrt{\frac{\tau_{w}}{\rho}}},$  where: ρ is the density ofthe fluid, and τ_(w) is the shear stress at the wall, expressed as${{\tau_{w} = {\rho\quad v\frac{\partial U}{\partial y}}}}_{y = 0}$  sothat the pressure velocity correlation can be expressed as:${- \frac{{{- \overset{\_}{p^{\prime}v^{\prime}}}}_{y = 0}}{{\rho U}_{\infty}^{3}}} = {{- \frac{1}{4}}\frac{H}{\delta}C_{k3}{\overset{\_}{u}}_{*}^{4}}$where:C_(k3)=C_(k2)K_(p)γ(ω) and where K_(p) is the Kraichnan parameter, witha value of 2.5, and γ(ω) is a dissipative function of the coatingmaterial, which reflects the phase shift between the pressurefluctuations and the vertical displacement of the coating, (c) the fluxof energy into the coating can be approximated by an effective turbulentdiffusion term, ${ɛ_{t}\frac{\partial k}{\partial y}},$  where:${ɛ_{t} = {C_{t}\frac{k}{ɛ}\overset{\_}{v^{\prime 2}}}},$  whereC_(t)=0.12  so that, in terms of nondimensionalized quantities:${{{- \frac{{{- \overset{\_}{p^{\prime}v^{\prime}}}}_{y = 0}}{{\rho U}_{\infty}^{3}}} \approx {{\overset{\sim}{ɛ}}_{t}\frac{\partial\overset{\sim}{k}}{\partial\overset{\sim}{y}}}}}_{y = 0}$ where: $\begin{matrix}{{\overset{\sim}{ɛ}}_{t} = \frac{ɛ_{t}}{U_{\infty}\delta}} \\{\overset{\sim}{k} = \frac{k}{U_{\infty}^{2}}}\end{matrix}$ and: $\overset{\sim}{y} = \frac{y}{\delta}$ (d) given theexpression for nondimensionalized pressure-velocity correlation in step(b) and in step (c), ε_(q)=ε_(t)|_(y=0) can be written as:${\overset{\sim}{ɛ}}_{q} = {\frac{ɛ_{t}❘_{y = 0}}{U_{\infty}\delta} = {\frac{1}{4}{C_{k3}\left( \frac{H}{\delta} \right)}{\overset{\sim}{u}}_{*}\frac{{\overset{\sim}{u}}_{*}^{3}\Delta_{\max}}{k_{\max}^{+} - k_{q}^{+}}}}$ where γ_(max) ⁺ is the normal distance from the surface to the maximumof turbulence energy, nondimensionalized as follows: $\begin{matrix}{y_{\max}^{+} = {\frac{y_{\max}}{v}u_{*}}} \\{\Delta_{\max} = \frac{y^{+}}{{Re}_{*}}}\end{matrix}$ where: ${Re}_{*} = \frac{u_{*}\delta}{v}$  is the Reynoldsnumber based on friction velocity, u*, and boundary layer thickness, δu* is the friction velocity, defined as μ₀ , and boundary layerthickness, δ μ₀ is the friction velocity, defined above as${u_{*} = \sqrt{\frac{\tau_{w}}{\rho}}},$ where: ρ is the density of thefluid, and τ_(w) is the shear stress at the wall, expressed as${{\tau_{w} = {\rho\quad v\frac{\partial U}{\partial y}}}}_{y = 0}$ and$k_{\max}^{+} = \frac{k_{\max}}{U_{\infty}^{2}}$  is the maximum ofturbulence kinetic energy, non-dimensionalized by U_(∞) ²$k_{q}^{+} = \frac{k_{q}}{U_{\infty}^{2}}$  is the kinetic energy of theoscillating surface, non-dimensionalized by U_(∞) ² (e) the boundarycondition for isotropic dissipation rate is expressed dimensionally as:${{{ɛ}_{y = 0} = {\frac{\partial}{\partial y}\left\lbrack {\left( {v + ɛ_{t}} \right)\frac{\partial k}{\partial y}} \right\rbrack}}}_{y = 0}$ or, in nondimensional form:${{{\overset{\sim}{ɛ}}_{y = 0} = {\frac{\partial}{\partial y}\left\lbrack {\left( {\frac{1}{Re} + {\overset{\sim}{ɛ}}_{t}} \right)\frac{\partial\overset{\sim}{k}}{\partial\overset{\sim}{y}}} \right\rbrack}}}_{y = 0}$where: ${Re} = \frac{U_{\infty}\delta}{v}$  is the Reynolds number basedon freestream velocity and boundary layer thickness.
 8. The method ofclaim 7, where the solution of the turbulent boundary layer problem foran isotropic viscoelastic surface is based upon an iterative technique,and where initial values are assumed for u* and for the gradient ofturbulent kinetic energy, ${\frac{\partial k}{\partial y}}_{y = 0}$ asbased upon the solution of the turbulent boundary layer problem for arigid flat plate, with two-dimensional flow, where the boundaryconditions are: $\begin{matrix}{{\frac{\overset{\_}{u^{\prime 2}}}{U_{\infty}^{2}}}_{y = 0} = 0} \\{{\frac{\overset{\_}{v^{\prime 2}}}{U_{\infty}^{2}}}_{y = 0} = 0} \\{{\frac{\overset{\_}{w^{\prime 2}}}{U_{\infty}^{2}}}_{y = 0} = 0} \\{{\frac{\overset{\_}{u^{\prime}v^{\prime}}}{U_{\infty}^{2}}}_{y = 0} = 0} \\{{{ɛ}_{y = 0} = {v\frac{\partial^{2}k}{\partial y^{2}}}}}_{y = 0}\end{matrix}$ where the values of u* and the gradient of turbulentkineticwhere the values of μ₀ and the gradient of turbulent kineticenergy, ${\frac{\partial k}{\partial y}}_{y = 0},$ obtained from thissolution are used to determine boundary conditions for the nextiteration, and where this procedure continues until the solutionconverges.
 9. The method for estimating the reduction in friction dragaccording to claim 8, for the case of two-dimensional flow over a flatplate coated by an isotropic viscoelastic material, wherein: a.) thevalue of the shear modulus, μ(ω), is constant in all directions, b.) thephase difference between normal and longitudinal displacements at thewall is equal to π/2, so that the Reynolds shear stresses at the surfaceof the viscoelastic coating are equal to zero:− u′v′=0 c) Reynolds stress boundary conditions are expressed innon-dimensional form as: $\begin{matrix}{\frac{\overset{\_}{u^{\prime 2}}}{U_{\infty}^{2}} = {\frac{1}{2}\left( \frac{H}{\delta} \right)^{2}\left( {C_{k1} + {{\overset{\sim}{u}}_{*}{Re}_{*}C_{k2}}} \right)^{2}{\overset{\sim}{u}}_{*}^{4}}} \\{\frac{\overset{\_}{v^{\prime 2}}}{U_{\infty}^{2}} = {\frac{1}{2}\left( \frac{H}{\delta} \right)^{2}C_{k2}^{2}{\overset{\sim}{u}}_{*}^{4}}} \\{\frac{\overset{\_}{w^{\prime 2}}}{U_{\infty}^{2}} = {\frac{1}{2}\left( \frac{H}{\delta} \right)^{2}C_{k1}^{2}{\overset{\sim}{u}}_{*}^{4}\tan^{2}\Theta}} \\{{- \frac{\overset{\_}{u^{\prime}v^{\prime}}}{U_{\infty}^{2}}} = 0}\end{matrix}$ where ${Re}_{*} = \frac{u_{*}\delta}{v}$  is the Reynoldsnumber based on friction velocity, u*, and boundary layer thickness, δu* is the friction velocity, defined as μ₀ , and boundary layerthickness, δ μ ₀ is the friction velocity, defined above as${u_{*} = \sqrt{\frac{\tau_{w}}{\rho}}},$  where: ρ is the density ofthe fluid, and τ_(ω)is the shear stress at the wall, expressed as$\begin{matrix}\begin{matrix}{{\tau_{w} = {\rho\quad v\frac{\partial U}{\partial y}}}}_{y = 0} \\{{\overset{\sim}{u}}_{*} = \frac{u_{*}}{U_{\infty}}} \\{C_{ki} = \frac{{\rho\quad K_{p}u_{*}^{2}{\xi_{i}}}❘_{p_{rms} = 1}}{H}}\end{matrix} \\{K_{p} = 2.5}\end{matrix}$ and where C_(k1) and C_(k2) can be approximated as zero ifthe amplitude of surface oscillations is less than the thickness of theviscous sublayer.
 10. The method according to claim 8, wherein thecoating is composed of multiple layers of isotropic viscoelasticmaterials with different material properties, where: a.) the shearmodulus, μ(ω), is constant within each of the multiple layers; b.)boundaries between layers are fixed, and c.) the static shear modulus ofthe material is progressively lower for each layer from the top layer tothe bottom layer of the coating.
 11. The method according to claim 8,wherein the coating is composed of an anisotropic material, where theproperties of the material in the normal direction (y) differ from thosein the transverse (x-z) plane, wherein the shear modulus in thestreamwise and transverse directions is given by μ₁(ω), and the shearmodulus in the normal direction is given by μ₂(ω): $\begin{matrix}{{\mu_{1}(\omega)} = {\mu_{01} + {\mu_{s1}\left\lbrack {\frac{\left( {\omega\tau}_{s1} \right)^{2}}{1 + \left( {\omega\tau}_{s1} \right)^{2}} + {{\mathbb{i}}\frac{{\omega\tau}_{s1}}{1 + \left( {\omega\tau}_{s1} \right)^{2}}}} \right\rbrack}}} \\{{\mu_{2}(\omega)} = {\mu_{02} + {{\mu_{s2}\left\lbrack {\frac{\left( {\omega\tau}_{s2} \right)^{2}}{1 + \left( {\omega\tau}_{s2} \right)^{2}} + {{\mathbb{i}}\frac{{\omega\tau}_{s2}}{1 + \left( {\omega\tau}_{s2} \right)^{2}}}} \right\rbrack}.}}}\end{matrix}$
 12. The method of claim 1, wherein step II thereofincludes the following substeps, not necessarily performed in the orderindicated: (a) choosing a density, ρ_(s), for the coating which iswithin 10% of the density of water, (b) selecting an initial choice forthe static shear modulus, μ₀, of the material to avoid resonanceconditions, based on the criterion that the speed of shear waves in thematerial, $\sqrt{\frac{\mu_{0}}{\rho_{s}}},$  is approximately the sameas the phase speed, C, of energy-carrying disturbances, (c) selecting aninitial choice for coating thickness, H, (d) solving the momentumequation for a viscoelastic material:${\rho_{s}\frac{\partial^{2}\xi_{i}}{\partial t^{2}}} = \frac{\partial\sigma_{ij}}{\partial x_{j}}$ given harmonic loading of unit amplitude, phase velocity, C,corresponding to the load of the turbulent boundary layer, and variablewavenumber and for a coating which is fixed at its base to a rigidsubstrate, where: ξ₁ and ξ₂ are the longitudinal and normaldisplacements through the thickness of the coating, approximated by thefirst mode of a Fourier series as follows:ξ_(i)=a_(i1)e^(ia) ^(e) ^((x−Ct)) where a_(i1) is a coefficient,α_(e)=ω_(e)/C is the wavenumber corresponding to maximum energy in theboundary layer, and φ_(i) is a phase difference, σ_(ij) is the stresstensor, written for a Kelvin-Voigt type viscoelastic material as:σ_(ij)=λ(ω)ε^(s)δ_(ij)+2μ(ω)ε_(ij) ^(s)  where ε_(ij) ^(s) is the straintensor, expressed as:$ɛ_{ij}^{s} = {\frac{1}{2}\left( {\frac{\partial\xi_{i}}{\partial x_{j}} + \frac{\partial\xi_{j}}{\partial x_{i}}} \right)}$ where ε^(s)=ε_(ii) ^(s), λ(ω) is the frequency-dependent Lame constant,defined in terms of the bulk modulus, K(ω), and the complex shearmodulus, λ(ω):${\lambda(\omega)} = {K_{0} - {\frac{2}{3}{\mu(\omega)}}}$ μ(ω) isconstant in all directions for an isotropic material, but will havedifferent values in different directions for an anisotropic material, and where the momentum equation is solved for a series of materialswhose shear moduli can be approximated for a single-relaxation type(SRT) of material, characterized mathematically in the form:${\mu(\omega)} = {\mu_{0} + {\mu_{s}\left\lbrack {\frac{\left( {\omega\tau}_{s} \right)^{2}}{1 + \left( {\omega\tau}_{s} \right)^{2}} + {{\mathbb{i}}\frac{{\omega\tau}_{s}}{1 + \left( {\omega\tau}_{s} \right)^{2}}}} \right\rbrack}}$ with an initial value of static shear modulus, μ₀=μ|_(ω=0), determinedin step 3(b), and for different values of dynamic shear modulus, μ_(s),and relaxation time, τ_(s), where:μ_(s)=μ|_(ω=∞)−μ₀ (e) from the series of SRT materials for whichcalculations were performed in step 3(d), selecting a material where thecalculated oscillation amplitude does not exceed the viscous sublayerthickness under given flow conditions, and where the energy flux intothe coating is considered over a range of frequencies approximately onedecade above and below the energy-carrying frequency, ω_(e), (f)iterating steps 3(d) and 3(e) with different values of thickness andstatic modulus, and determine the optimal combination of materialproperties for a single relaxation time (SRT) type of material, (g)using the results from step 3(f) as a guideline in the frequency rangeof interest, solving the momentum equation for viscoelastic materialswhose complex shear modulus may be characterized in the form of aHavriliak-Negami, or HN, material, whose shear modulus is expressed inthe following form:$\frac{\mu - \mu_{\infty}}{\mu_{0} - \mu_{\infty}} = \frac{1}{\left\lbrack {1 + \left( {\mathbb{i}\omega\tau}_{HN} \right)^{\alpha_{HN}}} \right\rbrack^{\beta_{HN}}}$where μ_(∞)=μ|_(ω=∞), τ_(HN) is a relaxation time, and α_(HN) and β_(HN)are constants, and (h) selecting final material properties of a coatingmaterial based on conditions of minimal oscillation amplitudes andmaximum energy flux, which is equivalent to the correlation betweenpressure and velocity fluctuations, − p′v′, in the range of frequenciescorresponding to maximum-energy-carrying disturbances, as well as basedon consideration of material fabricability.
 13. The method of claim 1,wherein step III thereof includes the following substeps, performed inthe order indicated: (a) solving the equations of motion and continuityfor a body with a viscoelastic coating, given the same location and flowconditions as in step I, using a numerical methodology which accountsfor non-zero energy flux and surface oscillation boundary conditions, aswell as for redistribution of energy between fluctuating components inthe near-wall region, and (b) comparing the friction drag calculated fora viscoelastic plate with that calculated for a rigid plate.
 14. Themethod of claim 1, wherein said viscoelastic coating is configured withmultiple surface or internal structures that approximate the naturallongitudinal and transverse wavelengths in the nearwall flow.
 15. Themethod of claim (1), wherein said viscoelastic coating is structuredwith a rigid underlying wedge or contour shape to minimize coatingthickness and oscillations near the intersection of the coating with arigid surface.
 16. A method for reducing the friction drag of a body byproviding a viscoelastic coating, said method comprising the followingsteps, performed in the order indicated: a) determining characteristicsof a turbulent boundary layer at a specified freestream velocity usingboundary conditions for a rigid surface of an identical size and shapeof the surface to be coated, said characteristics including viscoussublayer thickness, boundary layer thickness, phase speed and frequencycorresponding to maximum-energy-carrying disturbances, mean velocityprofiles, Reynolds stress distributions, wall shear stress distributionand friction drag; b) selecting material properties of the viscoelasticcoating including a density, a complex shear modulus and a coatingthickness that, when excited by a forcing function substantiallyidentical to that produced by the turbulent boundary layer as determinedin step a), will provide the maximum flux of energy into theviscoelastic coating without producing surface excitation amplitudesthat exceed the viscous sublayer thickness; c) determiningcharacteristics of the turbulent boundary layer over the viscoelasticcoating and at the specified freestream velocity using oscillationamplitudes and an energy flux corresponding to the material propertiesof the viscoelastic coating selected in step b), including mean velocityprofiles, Reynolds stress distributions, wall shear stress distribution,and friction drag, d) determining a percentage reduction in frictiondrag as a ratio of the value of the friction drag without theviscoelastic coating as determined in step a) minus the friction dragwith the viscoelastic coating as determined in step c) divided by thefriction drag without the viscoelastic coating as determined in step a),and thereby quantifying the performance potential of the coating design;and e) composing the viscoelastic coating from a material or acombination of materials as selected in all of the steps a) through d).17. The method according to claim 16, wherein the viscoelastic coatingis composed of multiple layers of isotropic, viscoelastic materials withdifferent material properties, where: a) the complex shear modulus,μ(ω), is constant within each of the multiple layers, b) boundariesbetween layers are fixed, and c) the static shear modulus μ ₀ of thematerial is progressively lower for each layer from the top layer to thebottom layer of the viscoelastic coating.
 18. The method according toclaim 16, wherein the coating is composed of an anisotropic material,and the properties of the anisotropic material in the normal direction(y) differ from those in the transverse (x-z) plane.
 19. The method ofclaim 16, wherein said viscoelastic coating is structured with a morerigid underlying wedge or a contour shape to decrease the viscoelasticcoating thickness and thereby minimize oscillations near theintersection of the viscoelastic coating with a rigid surface.
 20. Themethod of claim 16, wherein the viscoelastic coating is combined withsurface structures to enhance the stabilization of natural longitudinaland transverse wavelengths in the near-wall flow.
 21. A viscoelasticcoating designed according to the method of claim
 16. 22. A viscoelasticcoating designed according to the method of claim
 17. 23. A viscoelasticcoating designed according to the method of claim
 18. 24. A viscoelasticcoating designed according to the method of claim
 19. 25. A viscoelasticcoating designed according to the method of claim 20.